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The analogs of r-modes in superfluid neutron stars are studied here. These modes, which are 
governed primarily by the Coriolis force, are identical to their ordinary-fluid counterparts at the 
lowest order in the small angular-velocity expansion used here. The equations that determine 
the next order terms are derived and solved numerically for fairly realistic superfluid neutron-star 
models. The damping of these modes by superfluid "mutual friction" (which vanishes at the lowest 
order in this expansion) is found to have a characteristic time-scale of about 10 4 s for the m = 2 
r-mode in a "typical" superfluid neutron-star model. This time-scale is far too long to allow mutual 
friction to suppress the recently discovered gravitational radiation driven instability in the r-modes. 
However, the strength of the mutual friction damping depends very sensitively on the details of the 
neutron-star core superfluid. A small fraction of the presently acceptable range of superfluid models 
have characteristic mutual friction damping times that are short enough (i.e. shorter than about 
5 s) to suppress the gravitational radiation driven instability completely. 

PACS Numbers: 04.40. Dg, 97.60. Jd, 97.10. Sj, 04.30.Db 



I. INTRODUCTION 

Recently Andersson |ÏJ and Friedman and Morsink ||] 
showed that the r-modes in all rotating stars would be 
driven unstable by the emission of gravitational radiation 
in the absence of internal fluid dissipation. Subsequent 
analysis by Lindblom, Owen, and Morsink and then 
by Andersson, Kokkotas, and Schutz showed that in- 
ternal fluid dissipation in hot young neutron stars is in- 
sufficiënt to suppress this gravitational radiation driven 
instability. Thus neutron stars that are formed rapidly 
rotating are expected to spin down within about one year 
to a relatively small angular velocity (about 5-10% of 
the maximum) by the emission of gravitational radia- 
tion. Owen et al. || constructed rough models of this 
spindown process, and concluded that the gravitational 
radiation from these spindown events might be observ- 
able by the second-generation LIGO gravitational wave 
detectors. 

The purpose of this paper is to investigate the be- 
havior of this gravitational-radiation instability in the 
r-modes of older colder neutron stars Here the 

physics is more complicated and there are interesting ob- 
servational constraints. The existence of the two 1.6 ms 
pulsars and numerous examples of somewhat more 
slowly rotating neutron stars in low mass x-ray binaries 
(LMXBs) |Ï(J show that somc neutron stars are in fact 
rapidly rotating and stablc. Unfortunatcly uncertainty in 
the neutron-star equation of state means that the mini- 
mum rotation periods for neutron stars are not presently 
well known, and so we can not say exactly how rapid 
(in a dynamically meaningful sense) these rotations re- 
ally are. Values for the minimum rotation periods (when 
mass shedding first occurs) of 1.4 Mq models range from 



0.5 ms to about 1.4 ms, depending on the equation of 
state jnj. Any value in this range, however, implies that 
the 1.6 ms pulsars are rotating much more rapidly than is 
consistent with our present understanding of the gravita- 
tional wave spindown due to the r-mode instability in hot 
young neutron stars. It is widely believed, however, that 
these 1.6 ms objects are old cold recycled pulsars [ |Ï2| , 
having been spun up by accretion long after their ini- 
tial cool-down. These neutron stars are expected to have 
superfluid cores, and hence the fluid dynamics and dissi- 
pation mechanisms that govern their r-modes are entirely 
different from those studied to date The purpose 

of this paper is to develop the tools needed to study the 
r-modes in superfluid neutron stars. The challenge is to 
understand how the r-mode instability is suppressed in 
the 1.6 ms pulsars in particular, and the more numerous 
3 ms objects in LMXBs more generally. 

The superfluid dissipation mechanism called "mutual 
friction" seems a likely candidate to provide the needed 
stability for the r-modes in old cold neutron stars. Mu- 
tual friction arises from the scattering of electrons off 
the magnetic fields entrapped in the cores of the super- 
fluid neutron vortices |Ï1^ , |Ï4| ], and is known to play an 
important rolc in other aspects of the dynamics of su- 
perfluid neutron stars. Lindblom and Mendell |D| show, 
for example, that mutual friction completely suppresses 
the gravitational radiation driven instability in the ƒ- 
modes of rotating neutron stars. Our results here for the 
r-modes present a more ambiguous picture. We find in 
Sec. VI that the characteristic damping time for the r- 



modes due to mutual friction is about 10 s, for a typical 
model of the neutron-star core superfluid. This timescale 
is far too long to have any appreciable effect on the r- 
mode instability in these stars. However, we also find 
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that the mutual-friction damping time is extremely sen- 
sitive to the parameters that define the core superfluid. 
Within the presently acceptable range of the parameters 
examined here, about 1% have mutual friction damp- 
ing times so short (i.e. shorter than about 5 s) that 
the r-mode instability is suppressed completely. A some- 
what larger fraction of these parameters in the acceptable 
range, about 3%, have damping times short enough (i.e. 
shorter than about 58 s) that mutual friction suppresses 
the instability in some sufficiently warm and sufhciently 
slowly rotating neutron stars. Thus we conclude that 
an appropriately fine-tuned superfluid dynamics could 
provide the needed stability for the r-modes in old cold 
neutron stars through mutual friction. However given 
the small fraction of superfluid models that provide the 
needed stability, other damping mechanisms need to be 
considered (e.g. solid crust effects, strange quark matter, 
magnetic fields, ...). 

Regular shear viscosity is another mechanism that 
could play a role in damping the r-modes of superfluid 
neutron stars. We analyze this damping mechanism for 
the superfluid r-modes in Sec. |Vl| and find a characteristic 



damping time of about 10 (T/ÏO K) s, where T is the 
temperature of the neutron-star core. This time-scale is 
short enough to suppress the r-mode instability in stars 
cooler than about 10 6 K. Neutron stars are spun up in 
the usual picture during an LMXB phasc, in which the 
core temperature is expected to exceed 10 8 K jÏ6). This 
temperature is too hot to allow shear viscosity to provide 
the needed stability even for the 3 ms neutron stars ob- 
served in these systems. Levin || has shown that when a 
neutron star is spun up to the point where stability of the 
r-modes is lost, the star heats up and then spins down to 
a very small angular velocity in a few months by emit- 
ting gravitational radiation. Thus, some robust internal 
fluid dissipation mechanism must be identified to explain 
the stability of the r-modes in the observed LMXB sys- 
tems. If mutual friction is the only mechanism capable of 
providing the needed stability, then this fact would place 
interesting constraints on the parameters of the neutron- 
star core superfluid. An alternate possibility in the case 
of the 1.6 ms pulsars would be a mechanism for spinning 
up these stars without raising their core temperatures 
above about 10 7 K (e.g. by accretion at very low rates). 

In Sec. U we review the basic hydrodynamics of 
neutron-star core superfluids. We outlinc in Sec. III the 



derivation of the equations that govern the normal modes 
of a superfluid neutron star from this hydrodynamic the- 
ory. In Sec. |^ we take the small angular velocity ex- 
pansions of these equations that are needed to study the 
r-modes. In Sec. [v| we present our numerical solutions 
for the r-modes of rotating superfluid neutron stars, up 
to the second order in the small angular velocity expan- 
sion. The effects of superfluid mutual friction and shear 
viscosity on these r-modes are evaluated in Sec. VI. The 



II. SUPERFLUID HYDRODYNAMICS IN 
NEUTRON-STAR MATTER 

When the core temperature of a neutron star drops be- 
low about 10 9 K, a phase transition to a superconducting- 
superfluid state is expected to occur |Ï7|-|2"Ö[] . The neu- 
trons in the core are expected to form^-2 Cooper pairs 
and the protons to form 1 Sq pairs. The purpose of this 
scction is to review bricfly the equations that describe the 
behavior of this complicated superconducting-superfluid 
mixture on the macroscopic scales needed here to de- 
scribe the r-modes. 

Let v n denote the velocity of the neutron superfluid, 
and v p the velocity of the proton superfluid. On small 
scales these superfluid velocities, v n and v p , are related 
by the London equations to the phases of the complex 
order parameters, S n and S p , that describe the neutron 
and proton condensates; thus, v„ — (h/2m n )\? S n and 
Ü p = (h/2m p )V S p — (e/m p c)A, where A is the electro- 
magnetic vector potential. These equations imply that 
vorticity and magnetic fields in this material are confmed 
to vortices and flux tubes of microscopie dimension. In 
a typical neutron star the spacing between these neu- 
tron vortices is expected to be of order 10 -3 cm, while 
the spacing between magnetic flux-tubes is expected to 
be about 10~ 10 cm |ï^| . Our interest here is the very 
large-scale motions of this material associated with the 
low-order r-modes. Thus, it is appropriate to considcr 
all physical quantities, such as v n and v p , to be averaged 
over many vortices. The procedure for making this av- 
erage is described more fully by Bekarevich and Khalat- 
nikov Q, Baym and Chandler p||, Sonin [J||, Mendell 
and Lindblom Q , and Mendell |Ï5|,|6) . Throughout the 
remainder of this paper, all quantities are assumed to be 
so averaged. 

One of the interesting and unusual features of the 
neutro n-star co re superfluid is the so called "drag ef- 
fect" 
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equations that determine the superfluid pulsations are 
expressed in spherical coordinates in the Appendix. 



or "entrainment effect" |28[ . This ef- 
fect is caused by the fact that the conserved partiele 
currents are not simply proportional to the superfluid 
velocities. Instead these conserved currents are linear 
combinations of v n and v p : p nn Vn + Pn P v p for the neu- 
tron current, and p P pV p + p np v n for the proton current. 
Thus a given neutron superfluid flow v n is accompa- 
nied by a certain (small) current of protons, and vice 
versa. The mass-density matrix elements p nn , p pp , and 
p np are determined by the micro-physics of the many- 
body strong-interactions that occur between the neutrons 
and protons. This entrainment effect plays a crucial role 
in mutual friction (perhaps the most important dissipa- 
tion mechanism in this material) which we discuss in 
more detail in Sec. VI. Unfortunatcly these mass-matrix 
elements are not well determined at the present time. 
These quantities are constrained by Galilean invariance: 
Pnn = Pn~ Pn P and ppp — p p - p np , where p n and p p are 
the neutron and proton mass densities. But, the inde- 
pendent element p np must be determined directly from 
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the micro-physics. We find it convenient to re-express 
p np in terms of the dimensionless entrainment parameter 
e: 



fields in terms of Sv and Sw: the average and relative 
velocities of the core superfluids. These quantities are 
defined as 



Pnp 



-epn 



(2.1) 



Borumand, Joynt, and Kluzniak p9| estimate that e ~ 
0.04, and that its value is known at present only to within 
about a factor of two. Given this uncertainty we explore 
the properties of the r-modes over the expected range of 
superfiuid models with 0.02 < e < 0.06. 

The material in the core of a neutron star is a compli- 
cated mixture of neutrons, protons, electrons, muons, etc. 
While the general equations that describe the dynam- 
ics of this kind of charged superconducting-superfluid 
mixture have been studied |24[ , these general equations 
are considerably more complicated than are needed here. 
Our present interest is the dynamics of the superfiuid 
analogs of the r-modes; thus we are interested in dy- 
namics having length-scales comparable to the size, and 
time-scales comparable to the rotation period, of the 
star. Under these conditions the dynamics of the core su- 
perfiuid material simplifies considerably. On time-scales 
longer than the plasma time-scale (about 10~ 21 s) and 
the cyclotron time-scale (about 10~ 15 s), for examplc, 
this material is well described by the magnetohydrody- 
namic limit of the exact equations p6fl . In this limit 
the electrons and muons maintain exact charge neutral- 
ity with the protons. Similarly the electrons and muons 
are forced by scattering to move together as a single huid 
on time-scales longer than about 10 -9 s [^9|. Further, 
for dynamics on the time-scales of interest here, the bulk 
electrical currents are extremely small [^6|. Thus it is 
appropriate to simplify further and require the charged 
species to move together without generating any elec- 
trical current. The dynamical degrees of freedom of this 
material are reduced therefore to a pair of velocity vector 
fields — one for the neutrons and one for the protons — and 
a corresponding pair of thermodynamic scalar densities. 

In general there are forces in the complete dynamical 
equations (even for this reduced system) that describe 
the interactions between the smoothed superfiuid flow 
and the sheaf of vortices These additional forces 

are negligible for fluid motions with time-scales compa- 
rable to the r-modes ||^J2^], and we neglect them here. 
And finally, in general the dynamics would also include 
a "normal" component of this superfiuid material; how- 
ever, again we simplify by assuming that the temperature 
is well below the superfiuid transition and ignore these 
additional dynamical degrees of freedom. 

Our study is directed toward an exploration of the su- 
perfiuid analogs of the r-modes. We are primarily inter- 
ested therefore in examining the equations that describe 
the evolution of small departures from a uniformly ro- 
tating equilibrium neutron star. The dynamics of the 
neutron-star core superfiuid is described by two veloc- 
ity vectors and two thermodynamic scalar fields. It will 
be convenient to express the equations for these velocity 



and 



pSv = p n Sv n + PpSHp, 



Sw = Silp — Sv n . 



(2.2) 



(2.3) 



We use the prefix 8 to denote a small (Eulerian) per- 
turbation away from the equilibrium value of a quantity; 
while, quantities without prefix, such as p — p n + p p , 
denote the equilibrium values. The superfiuid vel ocity 
fields Sv n and Sv p are easily determined from Eqs. (2.2) 
and (2.3) once Sv and Sw are known. Similarly the ve- 
locity field of the electrons, 8v e , can can be expressed in 
terms of these quantities from the condition that there is 
no electrical current |2ql: 



Sv e = ^^-Sv p + -^-öv n = Sv + — 



Pp 



P 



e— ) Sw. (2.4) 
Pp/ 



We note that for simplicity in these equations we have 
ignored terms of order m e /m p , the ratio of electron to 
proton mass. And for simplicity in this discussion, we 
have also ignored the presence of muons. A more com- 
plete discussion including the contributions of the muons 
is given by Mendell [ p5| . 

The equations for the velocity fields Sv and Sw are ob- 
tained by perturbing the full system of superconducting- 
superfluid evolution equations, subject to the assump- 
tions described above. These equations, when expressed 
in terms of Sv and Sw, have the remarkably simple 
forms M: 



d t Sv a + v b V b Sv a +Sv b V b v a 



1 (dp 
~i? , .. . . 



*0V°p, (2.5) 



d t Sw a + v b V b Sw a + (27 - l)Sw b V b v a = -V a Sf3. (2.6) 

The pe rtur bed scalar SU that appears on the right side 
of Eq. (|2.5|) is defined by 



SU 



Sp 



- <5$, 



(2.7) 



where Sp is the perturbed pressure, and (5$ the perturbed 
gravitational potential . Th e potcntial S[3 that appears on 
the right sides of Eqs. (2.5) and (2.6) measures the degree 
to which the perturbed fluid départs from /3-equilibrium. 
The thermodynamic function f3 is related to the chemical 
potentials (per unit mass) of the neutrons //„, protons p, p , 
and electrons p, e by 



(3 = fi p — fi n -\ p, e 



(2.8) 
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The quantity (3 vanishes in the equilibrium state. (For 
simplicity we again neglect terms of or der m e /m v .) Fi- 
nally, the velocity field v a in Eqs. ( |2.5| ) and ( |2.6| ) repre- 
sents the uniform rotation of the equilibrium sta r, a nd 
the dimensionless quantity 7 that appears in Eq. ( [2.C| ) is 
related to the determinant of the superfluid mass-dcnsity 
matrix: 7 = (p nn p P p ~ p 2 np )/p n p p . 

The three perturbed scalar fields , SU , S$, and 8(3, that 
appear on the right sides of Eqs. (2.5) and ( |2.6D can be 
used to determine all of the other scalars of interest in 
this problem; for examplc 



6p = p(g) (8U + 6*) + (!£) 8 P , 



(2.9) 



and Spp (and thence öp n as well) can be determined from 



6(3 = 



dpp 
dp n 



Pp 



(Sp - Sp p ) 



+ 



dpp 
dp p 



dp n \ ml dp e 

9 Pp ) p„ m p d p e 



Sp p . (2.10) 



It is straightforward then to transform the mass conserva- 
tion laws (for neutrons and protons) into forms depend- 
ing only on SU, 5$, and 8(3 M: 



d t Sp + v a V a 8p + V a (p8v a ) = 0, 



dp 

djü 



(SU + 8$) 



p 8(3 



dp 

9(3 



8v a V aP + V a (p8w a )=0. 



(2.11) 



8(3 



(2.12) 



The quantity p that appears in Eq. ( 2.1 2f ) is defined as 
P = (P nnPvv-pl J /p = PnPpl/p- These Eqs. ( pl)| ), (pTë|), 
(2.11), and ( [2.12 ), together with the perturbed gravita- 
tional potential equation, 



V a V Q 5$ = -AwGSp. 



(2.13) 



determine the evolution of the material in the superfluid 
core of a neutron star in the long length-scale, long time- 
scale, and low temperature approximation of interest to 
us here. 



III. OSCILLATIONS OF SUPERFLUID 
NEUTRON STARS 

A superfluid neutron star is a reasonably complicated 
structure consisting of a superfluid core surrounded by 
a solid crust, and probably a liquid ocean above that. 
For the purposes of our analysis here we use a simplified 
and idealized representation of this structure. We con- 
sider a neutron-star model that consists of a superfluid 
core (where p > p s ), surrounded by an ordinary matter 



envelope (where p < p s ). For simplicity we treat the ma- 
terial in this envelope as a perfect Huid. The dynamics of 
t he m at erial in the core is described by Eqs. (2.5), (|2.6|), 
( |2.1l|) , ( 2.12 ), and ( 2.13 ). And similarly the material in 
our idea lized envelope is descri bed b y Eul er's e quation, 
Eq. (|J) with 8(3 = 0, and Eqs. fë.ll| ) and ( fÜj ). In this 
section we show how the modes of such a rotating super- 
fluid stellar model can be described completely in terms 
of the three scalar potentials SU, 8(3 and 8$ Q. And 
we derive the equations that determine these potentials, 
together with the appropriate boundary conditions. 

We assume here that the time dependence of the per- 
turbation is e lwt and that its azimuthal angular depen- 
dence is é lmv , where ui is the frequency of the mode and 
m is an integer; The su perfl uid versions of the Euler 



equation, Eqs. (2.5) and (2.6), determine the velocities 
8v a and 8w a in terms of the scalars SU and 8(3 much as 
they do in the ordinary- fluid case pï| . Given the tem- 
poral and angular dependence assumed here, Eqs. (2.5) 
and (|j]) become linear algebraic equations for 8v a and 
8w a which can be solved directly: 



Sv a = iC 



6w a = iQ ab S7 b S(3. 



(3.1) 



(3.2) 



In these equations Q ab and Q ab are tensors that depend 
on the frequency of the mode u>, and the angular velocity 
of the equilibrium star fi. These tensors are given by: 



Q 



ab 



1 



(lü + mVL) 2 - 4Ü 2 

(üj + m n)s ab - 



Aü 2 z a z b 
lü + mfl 



2iOVV 



(3.3) 



Q 



ab _ 



(lü + mü) 2 - 4 7 2 f7 2 

4-/ 2 n 2 z a z b 

lü + mfl 



(lü + mQ)S ab 



2i7fïVV 



(3.4) 



In Eqs. ( p.3[ ) and (3.4) Q is the angular velocity of the 
equilibrium star; the unit vector z a points along the ro- 
tation axis; ip a is the vector field that generates rotations 
about the z a axis; and 8 ab is the Euclidean metric tensor 
(the identity matrix in Cartesian coordinates). 

The expressions for the velocity fields in Eqs. ( ft.l| ) and 
T2[) can be substituted into the mass conservation laws, 



Eqs. (2.11) and (2.12), to obtain equations for the scalar 
fields alonc. In general, the potentials 8U, 8(3 and <5$ are 
solutions then of the following system of partial differen- 
tial equations |3G]: 
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V a ( P Q ab W b SU) + (uj + mÜ)p ^ J SU 



Mm), 11 "** 



— (w + miï) 



Ö/3 



(3.5) 



V a (pQ o6 V 6 ö/3) + (w + mfi)^— 5/3 

-AfëT) Q afc v a pV oP 5/3 = 



ö/3 



9p 
95 



-Q ab V a pV b SU +(lü + mtt) (SU + 5$) 
P 



(3.6) 



V a V a< 5$ + 47rGp ( |h <5$ 

\ d Pj i3 



dp 



The functions SU, 5/3, and 5$ are also subject to ap- 
propriate boundary conditions at the interface between 
the superfluid core and the ordinary-fluid envelope, at the 
surface of the star, and at infmity. First, we consider the 
boundary at the interface between the superfluid core and 
the ordinary-fluid envelope of the star. Mass and momen- 
tum conservation across this boundary place a number of 
constraints on the continuity of these functions |3Ö|. In 
particular these conditions require that the functions SU 
and 5$ be continuous there. In addition, V a 5$ must be 
continuous, while V a SU must have a discontinuity that 
is prescribed by 

1 (dp 



d/3 



a [V a SU} s - _ n a V aP [Sf3} s = n a [V a 5U] . 



(3.8) 



The subscripts s and o in Eq. (3.8) denote that the quan- 
tities are to be evaluated as limits from the superfluid 
or the ordinary-fluid side of the boundary respectively, 
and n a denotes the outward directed unit normal to the 
boundary surface. The function 5/3, which is of interest 
to us only within the superfluid core, must satisfy the 
condition 



n a [V 5/3], 



- — ■ —-^z V a 5/3 S 

[ui + muy 

2m^Vlvj b n b 
(lu + mQ)zu 



0, 



(3.9) 



on the boundary of the superfluid core. Here we use the 
notation w for the cylindrical radial coordinate, and w a 
to denote the unit vector in the w direction. 



Next consider the boundary conditions on the outer 
surface of the star. The function SU must be constrained 
at this surface in such a way that the Lagrangian per- 
turbation in the pressure vanishes there: Ap = 0. This 
condition can be written in terms of the variables used 
here by noting that 



Ap = Sp ■ 



5v a ^aP 
i(u> + mfl) ' 



(3.10) 



where Sv a is given in this re gion by Eq. fl3.l| ) with 5/3 = 0. 
Thus using Eqs. (2.7) and (3.1) this boundary condition 
can be written in terms of SU and 5$ as 







p (u + mü)(SU + 5$) + Q ab V aP V b SU 



(3.11) 



Finally, the perturbed gravitational potential 5$ must 
fall off at infmity faster than l/r in order that the mass 
of the perturbed star remain the same as that of the equi- 
librium star: lim r _>oo (r5$) = 0. In addition 5$ and its 
first derivative must be continuous at the surface of the 
star. The problem of finding the modes of "uniformly" 
rotating superfluid st ars is red uced t here fore to finding 
the solutions to Eqs. (3.5), (3.É) and (3.7) subject to the 
appro pria te bou ndary cond itions including in particular 

Eqs. (|j), f[|), and CU)- 

The cquations for the potentials SU and 5/3, Eqs. (3.5) 



and (3.6), have complicated dependences on the fre- 
quency of the mode and the angular velocity of the star 
through Q ab and Q ab , as given in Eqs. and Q). In 
the analysis that follows it will be necessary to have thosc 
dependences displayed more explicitly. Here we are inter- 
ested in investigating the superfluid versions of r— modes. 
Such modes have frequencies that go to zero linearly as 
the angular velocity of the star vanishes. Thus, it will be 
uscful to define the dimcnsionless frequency parameter k: 



kU, — lu + mQ. 



(3.12) 



The parameter k remains finite in the zero angular- 
velocity limit for these modes. Using this para mete r an d 
the expressions for Q ab and Q ab from Eqs. Q) and 
we re-write Eqs. (3.5) and ( |3.6| ) to obtain the following 
equivalent forms: 



p{n 2 5 ab - Az a z b )V b 5U 
-k 2 (k 2 -4)ü 2 
1 (dp 



2mn 



w 



w a \7 aP SU = 



-v„ 

2rriK 



p \df3 



dp\ 



pzu \df3J p 



5/3 (K 2 S ab - 4z a z b )V bP 



Spm a V a p, 



(3.13) 



5 



+ 



K 2 — A"f 2 

IrriK 



( K 2 5 ab - 4j 2 z a z b )V b 6(3 



— 



1P 



K — 47 
2 



P df3 \p 



5(3 



dp 



K 2 S ab _ 4z a z b 



p(k 2 - 4) 



V a pV b 5U 



IrriK 



P(k 2 - A)vj 



TU a V a pSU + K 2 n 2 (SU + (5$) 



(3.14) 



The boundary conditions, Eqs. (3.9) and (3.11) are sim- 
ilarly transformed into the forms: 

K 2 n a [V a ö(3} s - ^ 2 z b n b z a ^ a 8(3] s + 2mKlwbnb [5(3] s = 0. 

zu 

(3.15) 



K 2 S ab -4z a z b )VahV b SU- 



2mn 



-w a V a h SU 



+k 2 {k 2 ~4)n 2 {5U + 6$) 



= 0. 



(3.16) 



In Eq. ( 3.16 ) we have expressed the boundary condition 
in terms of the thermodynamic enthalpy, h, which is de- 
fined as 



h(p) 



Wï 



(3.17) 



The enthalpy is the appropriate thermodynamic function 
to use in Eq. (3. IÉ) because its gradiënt, V a h, is well- 
defined and non-zero at the surface of the star. 



P = PO + P2 



nGpo 



h Q 



Ü 2 

> — — 

'nGpo 



(4.2) 



(4.3) 



The location of the surface of the star R(n) is also ex- 
pressed as such an expansion: 



R = i?o + R2 



9 2 



nGpo 



(4.4) 



Here and throughout the remainder of this paper we use 
the subscripts and 2 to denote the lowest- and the 
second-order terms respectively in these expansions; and 
we use r and \x = cos 9 to denote the Standard spherical 
coordinates. We also introducé here the angular veloc- 
ity scale \/nGpo, where po is the average density of the 
star in the non-rotating limit. (Equilibrium neutron-star 
models do not exist for f2 > |y / ^G/5Ö.) The techniques 
needed to evaluate the terms in these series for the equi- 
librium structure are identical to those described for ex- 
ample in Lindblom, Mendell, and Owen |33| . 

Next, we define expansions for the quantities that de- 
termine the perturbations of a superfluid star, SU, 5(3, 
(5$ and k: 



5U 



Rfo 2 



8(3 = Rltt 2 



5<I> = Rfo 2 



5U + 5U 2 



6 Po + 6fo 



n 2 



nGpo 



9 2 



nGpo 



(5$ + <5 < E>2 



a 2 



-n-Gpo 



ö{ü 4 



Ö{Ü A ) 



0{Q 4 



(4.5) 



(4.6) 



(4.7) 



IV. SUPERFLUID R MODES IN THE SLOW 
ROTATION APPROXIMATION 

The r-modes of rotating ordinary-fluid stars have tra- 
ditionally been studied using a small angular- velocity ex- 
pansion Our goal here is to perform a similar ex- 
pansion for the superfluid generalizations of th e r-m odes . 
Thus, we seek solutions to Eqs. (|j|), ( |3.14| ) and Q) 
as power series in the angular velocity of the star. 

We begin first with the structure of the equilibrium su- 
perfluid star. This structure is identical to its ordinary- 
fluid counterpart in the large-scale averaged-over-vortices 
hydrodynamics used here. We expand each of the equi- 
librium functions of interest: 



P = Po + P2 



tl 2 
irGpo 



(4.1) 



n 2 



K = K + K 2 



nGpo 



0{Q 4 



(4.8) 



We have normalized the eigenfunctions using Rq, the ra- 
dius of the star (in the non-rotating limit), and fi, the 
angular velocity of the star. Using these expressions for 
the perturbations, together with those for the structure of 
the equilibrium star, it is straightforward to write down 
order-by-order the equations that determine the super- 
fluid r-modes. This expansion is completely analogous 
to that given by Lindblom, Mendell, and Owen |33| for 
the ordinary-fluid modes. 

It is straightforward to verify that the functions 



SU 



r 

Rq 



ö(3 = 0, 



(4.9) 



(4.10) 
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with 



1 



(4.11) 



satisfy the lowe st-or der ter ms f rom the expansion of the 
pulsation Eqs. ( Ë>.13|) and ( |3.14 ). Inspection of the equa- 
tion for 5 f3, Eg. ( ft.!4| ), rcveals that the right side is pro- 
portional to f2 3 : the first ter m on the right vanishes iden- 
tically for 5Uq given by Eq. ( |4.9| ), while the lowest-order 
contribution to the second term is proportional to il 3 . 
Thus the function S(3q = satisfies Eq. ( 3.14 ) to lowest 
order. 

With the lowest-order contribution to 8(3 vanishing, 
the lowest-order equations for SU and <5$ from Eqs. ( 3.13j ) 
and (3/7) reduce to the following, 



p (i4ö ab - Az a z b )V b 5U Q 



— 



vj a V a po SU = 0, 



V a V a (5$ = -4ttG 



(SU + <5$ ) 



(4.12) 



(4.13) 



The lowest-order boundary conditions a t th e superfluid 
ordinary-fluid interface, Eqs. (3.5) and (3.9), merely re- 
quire that SUq and V a (5£/o are continuous there. The 
lowest-order contribution from the boundary condition 
for öUq at the surface of the star is 



{45 



-4z a z b )V a h V b SU Q 

IrriKn 
+ -vo a V a h ÓU 



= 0. 



(4.14) 



r=R 



These equations are identical the lowest-order terms in 
the ordinary-fluid r-mode equations [B3[. The function 



6Uq t ogeth er wit h kq given in Eqs. (4.9) and (4.11) satisfy 
Eqs. (4.12) and (4.14) identically because these are in fact 
the lowest-order expressions for the classical r— modes 
as studied for example by Papaloizou and Pringle [^2[ . 
(However, they are expressed here in a form that was 
introduced more recently j34j.) Thus to lowest order the 
superfluid r-modes are identical to their ordinary-fluid 
counterparts. 

Continuing on to second order, the equations for the 
potentials are 



V a 



p (4S ab ~ 4z a z°)V^C/ 2 



a b\ 



+ ^^v aP0 su 2 



p 2 (4ö ab - 4z a z")V b SU a + 2 K0 K 2 p a V a SU 



H Zü a (K2^aP0 + K V a p 2 ) SU 

VJ 

dp 



-TrGp 4(4 - 4)p TT {$Uq + M> ) 



dp 







Ml) ^ 



4z a z b )V bPo 



2m«o 1 f dp 



vj po \d/3 



5/3 2 w a \7 a po, 



(4.15) 



Po 



«o- 4 7o 



VJ a V a 



(4ö ab - ^z a z°)V b öp 2 
loPo 



4 ~ 475 

2 



Sf3 2 



pg(«§-4) \df3 

dp\ 

po(«§-4) \dll) 



% ) (4ö ab - 4z a z b )X7 a p X7 bPo Sfc 



(4ö ab - 4z a z b )V a p V b SU 2 



p 

^^-vj a \7 aPo 5U 2 + 2n K 2 X7 a p V a 5Uo 
vj 

-^l^V aP0 SU Q + ^vj«V aP2 SU 

VJ VJ 



+ (45 ab - 4z a z b )V a p 2 V b 5U 
-4irGp (^] (6U + ó$o) 



d/3 



(4.16) 



V"V>Ï>2 + 4^,, ( -J-j S<ï> 2 = 



dp\ 



-AirGpo 



Öp\ 
dp) . 



SU 2 - 4ttG 



dp\ 

dl3J v 



S(3 2 



dp) p dp 



Po 



dp\ 
dp) 



p 2 (óU + ö<ï>o). 



J /3 



(4.17) 



It will be helpful to express Eqs. ( 4.15| ) and ( 4.16| ) in 
the following shorthand forms: 



D{SU 2 ) =E{8p 2 )+F, 
D{6/3 2 ) = É(ÖU 2 ) + F, 



(4.18) 



(4.19) 



where D and D are second-order partial differential oper- 
ators, E and E are first-order operators, and F and F are 
functions that depend on the lowest-order cigenfunctions 
SU and <5$ . 

We next consider the second-order boundary condi- 
tions. The boundary conditions that must be satisfied 
at th e superfluid ordinary-fluid boundary, Eqs. (3.S) and 
(|3.9|) , have the following second-order forms on the sur- 
face r = R s : 



döU 2 



1 / dp\ dp a 
72 ( ) — [<%L = 



Po \dPJp dr 



dSU 2 



(4.20) 



(4 - 4 7o V) 



d5/3 2 
dr 
2mK 7o 



P{l-P 2 ) 



d8f3 2 
dp 



[Sfo], = 0. 



(4.21) 
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The seco nd-or der boundary condition for the potential 
SU, Eq. ( 3.16 ), is identical to that derived by Lindblom, 
Mendell, and Owen |33| for the perfect-fluid case: 

| {n 2 5 ab - Az a z b )V a h V b SU 2 + ^^w a V a h 5U 2 

+ { K 2 S ab - 4z a z b )V a h 2 V b ÖU + ^^nj a W a h 2 SU 
+2K K 2 S7 a h \7aöUo + ^— ^ a V a h 6U 



+k 2 (k 2 - 4)irGp Q (6U + 6<ï> ) 



(4.22) 



It will be useful in our numerical solution of the pulsa- 
tion equations to have an expression for the second-order 
contribution to the frequency of the mode, n 2 , as inte- 
grals over the eigenfunctions. Such an expression can be 
obtained by multiplying Eq. J4.1S| ) by SUq and integrat- 
ing over the interior of the star. Since SUq is an element 
of the kernei of the operator D, this part of the integral 
reduces to a boundary integral. This boundary integral 
is non-vanishing because V a SU 2 is discontinuous at the 
superfLuid bou ndary as a consequence of the boundary 
condition Eq. ( 4.2C| ). The result of integrating the left 



side of Eq. (4.18) can be expressed therefore in the fol- 
lowing way: 



SU;D{SU 2 )d 3 x 



2tt 



r=R B 

(4.23) 



Combining this result with the more straightforward in- 
tegrals of the right side of Eq. ( 4.1S ) , we obtain the fol- 
lowing expression for the second-order change in the fre- 
quency of the superffuid r-mode: 



«2 



--j-lëUo] d a X = —2 

r dr (m + l) z 



P22 



\6U \*d*x 



SirGpom 



dp 



+ 



(m + l) 4 

(k% - 4p 2 )r d 
2(m + 2) dr 



(SU + ö$ )6Uö* d 3 x 



1 dh a 
r dr 



dp 

djï 



Sp 2 SU^d ó x 



dh a ( dp \ 
dr \dp). 



3kq — 2mKo — 4 
2(m + 2) 



2 dho fdp\ 
dr {dHj, 



nl - 4p 2 dS/3 2 
2(m + 2) dr 

2p(l - p 2 ) 85p 2 
(m -I- 2)r dp 



SUod 3 x 



1 dh ( dp\ XTT *,3 



4p< 



Sp 2 SUgdp 



(4.24) 



In Eq. (4.24) the integrations are to be performed over 
the interior of the superfiuid core, < r < R s , for those 
integrals involving ö/3 2 and throughout the interior of the 
star, < r < i?o> for those integrals that do not. Even 
though this expression for k 2 depends on the second- 
order eigenfunction S(3 2 it is nevertheless very useful in 
determining k 2 numerically. This is due to the fact that 
the integrals involving 5j3 2 are typically quite small. 



V. NUMERICAL SOLUTIONS FOR THE 
SUPERFLUID fl-MODES 

In order to solve the equations for the superfiuid r- 
modes we must adopt a specific model for the equilibrium 
structure of the neutron star, and explicit expressions for 
the various thermodynamic functions that appear in the 
equations. For the equilibrium structure of the neutron 
star, we use the simple polytropic equation of state: p = 
Kp 2 , with K chosen so that a 1.4M Q model has a radius 
of Rq — 12.533 km. We cho ose t his simple model because 
our method of solving Eq. ( 4.18 ) to determine SU 2 seems 
to be rather unstable for more realistic equations of state. 
We also use this simple equilibrium equation of state to 
evaluate the thermodynamic derivative (dp/dp)p — p/2p 
that appears in the pulsation equations. We adopt the 
fiducial value p s = 2.8 x 10 14 g/cm 3 for the superfiuid 
transition density. However, we examine the sensitivity 
of our results to variations in this parameter. 

In order to evaluate the other thermodynamic proper- 
ties of neutron star matter needed in the pulsation equa- 
tions, we use the recent semi-empirical equation of state, 
A18+<5u+UIX, of Akmal, Pandharipande, and Raven- 
hall |3^]. The derivative {dp/d(3) p is determined from 
this equation of state using the formula, Eq. (70), in 
Lindblom and Mendell 1 30 1 that expresses this deriva- 
tive in terms of the more easily evaluated derivatives of 
the chemical potentials. For the density range of primary 
interest to us (2.8 x 10 14 < p < 10 15 gm/cm 3 ), we find 
that 



dp 



-h-) w 1.4x10" 



1.1 x 10 



'22 , 



(5.1) 



where all quantities are expressed in cgs units. Similarly, 
we find that the proton fraction p p / p for this equation of 
state is given approximately by 



Pp 
P 



0.031 + 8.8 x 1Q- I7 p. 



(5.2) 



These expressions are accurate to within about 20% in 
the indicated density range. 

The thermodynamic quantities p and 7 can be ex- 
pressed in terms of the baryon densities p n and p p , 
and the mass-density matrix element p np quite generally 
using the constraints on th e m ass-density matrix from 
Galilian invariance and Eq. (2.1): 



r=R s 
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Pp 



pi 



pp 
p 



Pp 
p 



(5.3) 



(5.4) 



Thus using Eq. ( 5.2 ) for the proton fraction, we obtain 
simple smooth expressions for t he su perflui d thc rmody- 
namic functions needed in Eqs. (4.18) and (4.19). 

In order to evaluate the potentials SU2 and (5/32, we 
convert the operators in Eqs. (4.18) and ( 4.19| ) into ma- 
trix operators by making the usual discrete representa- 
tions of the various derivatives that appear. To this end 
we have given expressions for these operators in sphcri- 
cal coordinates in the Appendix. We convert the radial 
derivatives into discrete form using the Standard three- 
point differencing formulae. For the angular derivatives, 
we have written separate codes that use either the higher- 
order angular differencing formulae discussed in Ipser and 
Lindblom |H| or the Standard three-point formulae. We 
find that the results of our two separate codes agree very 
well. The code that uses the higher-order angular dif- 
ferencing formulae requires far fewer angular spokes in 
order to achieve a given accuracy, and it appears to be 
somewhat more stable. 



The operator D that appears in Eq. (4.18) is a hyper- 
bolle differential operator with a non-trivial kernei. Thus, 
Eq. (4.18) is solved using the numerical relaxation tech- 
nique developed in Lindblom, Mendell, and O wen |53|. 
The operator D is also hyperbolic, but (generically) it 
does not have a non-trivial kernei. We find that D canbe 
inverted numerically without difnculty. Thus, Eq. ( 4.19| ) 
can be solved by straightforward numerical tcchniques. 

The Eqs. ( 1.18| ) and (4.19) are solved iteratively to 
determine the potentials 6U2 and 8(32- We begin the pro- 
cess by setting (5/32 = 0. Then we iterate the follow- 
ing sequence of calculations. First, we evaluate K2 from 
Eq. ( f4.24| ); t hen w e solve Eq. ( f4. 18] ) for 5U 2 ] and finally 
we solve Eq. ( 4.19 ) for 50 2 - This sequence of calculations 
is iterated until the functions ÓU2 and (5/32 converge. This 
generally takes only a few steps. We find that the second- 
order eigenvalue K2 has the value 0.29879 for this model, 
essentially independent of the entrainment parameter e 
and the superfluid-transition density p s . This value of «2 
differs only slightly from the value 0.29883 found by Lind- 
blom, Mendell, and Owen |53) for the analogous ordinary- 
fluid r-modes. 

The eigenfunction 5U2 for the case e = 0.04 and 
p s = 2.8 x 10 14 g/cm 3 is illustrated in Fig. [[]. This func- 
tion was computed on a grid having 2000 radial grid 
points and 10 angular spokes using the higher-order an- 
gular differencing code. Each curve represents the radial 
dependence of ÖU2 along one of the angular spokes. This 
function is essentially identical to that obtained by Lind- 
blom, Mendell, and Owen J33| for the ordinary-fluid case. 
Fig. U illustrates the associated function 502 that was 
evaluated together with the 6U2 of Fig. |ï|. We note that 



5U, 




FIG. 1. Eigenfunction 8U2 for the m — 2 r-mode based on 
a superfluid model with e = 0.04 and p s = 2.8 x 10 14 g/cm 3 . 
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8p 2 0.02 
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r/R 

FIG. 2. Eigenfunction 502 for the m = 2 r-mode based on 
a superfluid model with e = 0.04 and p a — 2.8 x 10 14 g/cm 3 . 



(5/32 is only defined within the core of the star where the 
neutron-star matter is superfluid. We also note that 502 
is about an order of magnitude smaller than 5U2- Thus 
even at the second-order, the superfluid r-modes diffcr 
little from thcir ordinary-fluid counterparts. 

Figures | and ^| illustrate how 502 changes as e varies. 
In particular Figs. ||| and |] illustrate 502 for the extreme 
values of e considered here: e = 0.02 and e = 0.06 re- 
spectively. Whilc the differences in these three func- 
tions, Figs. ^, H, and [|, are significant, they do not re- 
ally illustrate the most interesting feature of their de- 
pendence on e. The most interesting and unexpected 
feature that we find in the solutions for 502 is a kind of 
resonance phenomenon. We find that there are certain 
"critical" values of e (e c = 0.02294 and e c = 0.04817 for 
the p s = 2.8 x 10 14 g/cm 3 case) near which the function 
502 becomes extremely large. Near these special values 
of e the character of the second-order terms in the ex- 
pansion for the mode change from being dominated by 
the ordinary-fluid like correlated motion of the neutrons 
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FIG. 3. Eigenfunction S/32 for the m = 2 r-mode based on 
a superfluid model with e = 0.02 and p a = 2.8 x 10 14 g/cm 3 . 
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FIG. 4. Eigenfunction 5/32 for the m — 2 r-mode based on 
a superfluid model with e — 0.06 and p a = 2.8 x 10 14 g/cm 3 . 



and protons, to the uniquely superfluid anti-correlated 
motion, where the neutrons and protons have opposite 
velocities. We see a hint of this behavior, perhaps, in 
Fig. U where the magnitude of 6 fa for e = 0.02 is consid- 
erably larger than its value at e = 0.04 or 0.06. 

Mathematically we find that the reason for this res- 
onance phenomenon is that the operator D that deter- 



mines 8 fa in Eq. (4.19) becomes singular at these critical 
values of e. This singular behavior is illustrated in Fig. ||, 
in which we depict the dependence of the smallest (in ab- 
solute value) eigenvalue of the operator D as a function 
of e. We see that this smallest eigenvalue vanishes for 
the two critical values of e noted above: e c = 0.02294 
and e c = 0.04817. In a sufficiently small ncighborhood of 
these points the function 5 fa must be proportional there- 
fore to i/>min/A m i n where V'min is the eigenfunction of D 
that corresponds to this smallest eigenvalue, A min . Since 
Amin vanishes for the singular values of e, this forces the 
function 5 fa to become very large in sufficiently small 
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FIG. 5. Eigenvalue of the operator D having the smallest 
absolute value as a function of e, for m — 2 perturbations 
with p a = 2.8 x 10 14 g/cm 3 . 



neighborhoods of these points. While this vanishing of 
the eigenvalues of D explains mathematically why 5 fa 
becomes large for certain values of e, it does not really ex- 
plain physically the cause of this resonance phenomenon. 



VI. DISSIPATION IN SUPERFLUID iï-MODES 

The primary motivation for our study here is to in- 
vestigate how the transition to a superfluid state at low 
temperatures effects the stability of the r-modes in older 
colder neutron stars. The r-modes are driven towards 
instability by gravitational radiation but internal 

fluid dissipation tends to suppress this instability HQ. 
Internal fluid dissipation can stabilize the r-modes com- 
pletely, if it is sufficiently strong. In this section we inves- 
tigate the importance of the two types of internal fluid 
dissipation considered most likely to have a substantial 
influcncc on the r-modes in superfluid neutron stars: mu- 
tual friction caused by the scattering of electrons off the 
cores of the neutron vortices, and shear viscosity due to 
electron-electron scattering. 

The effects of dissipation on the evolution of a fluid 
are most conveniently studied using an appropriate en- 
ergy functional. For the case of a neutr on- st ar su p erflui d 
whose evolu tion i s determined by Eqs. (2.5), (2.6), ( 2.11 ), 
Q2.12 ), and (2T3), the following is the appropriate energy 
functional [B6l 



C 2 



\9pJp 

+pSw:Sw a + ^^-(^-) 6/3*6(3 

P OP \PnJ p 



dp 

d3 



Re[(2ÓU + 6^)5/3*} }d 3 x. (6.1) 
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The integrals of the first two terms on the right side of 
Eq. (3.1) are to be performed throughout the star, while 
the last three terms (those proportional to 6w a or 6/3) are 
to be integrated only within the superfluid core. In the 
absence of dissipation, the energy £ is conserved. When 
we evaluate t he s mall angular-velocity expansion of the 
energy in Eq. ( |6.l[ ) , we find that only the first term on the 
right contributes at the lowest order. Thus, the lowest- 
order expression for the energy, 



7 7r(m + l) 3 , 
£ = a 2 — — —^-(2m + 1)1 x 
2m 



fRo / \ 2m+2 



(6.2) 



is identical to its ordinary-fluid counterpart p3| . 

Mutual friction tends to damp out any relative motion 
between the proton and neutron superfluids. This dis- 
sipation is caused by the scattering of electrons (whose 
motion tracks that of the proton superfluid in our ap- 
proximation) off the magnetic fields that are entrapped 
within the cores of the neutron vortices. The presence 
of these magnetic fields is due to the entrainment effect, 
and hencc the magnitude of this dissipation mechanism 
depends strongly on the poorly known mass-matrix el- 
ement p np . The rate at which energy is dissipated by 
mutual friction in the neutron-star superfluid considered 
here is given by |3q| : 



d£_ 

dt 



2Ü / B nPn Y{5 a - z a z b )6w a 6w* b d i x, 



(6.3) 



where the dimensionless mutual friction scattering coëf- 
ficiënt, B n , is given by, 



2 1/6 

e Pp 



1.96 x 10 4 \p p 



-1 1-e 



Pp 



-3/2 



(6.4) 



We note that unlike regular viscosity, the mutual fric- 
tion scattering coëfficiënt does not depend on tempera- 
ture (at least at the lowest order). This independence is 
due to the fact that mutual friction involves the scatter- 
ing of two different fluids in relative motion. Thus, even 
at zero temperature some electrons acquire sufficiënt en- 
ergy from their collisions with the vortices to scatter into 
available energy states located above the Fermi surface. 

In order to evaluate the effects of mutual friction o n the 
superfluid r- modes, it is necessary to re-express Eq. (6.3) 



in terms of the eigenfunction 6/3 2 that determines the 
relative superfluid motion. To lowest order in the angular 
velocity then, we find 
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FIG. 6. Characteristic damping time t M f due to superfluid 
mutual friction. These curves show the dependence of f MF as 
a function of the superfluid entrainment parameter e for three 
values of the superfluid transition density p s - 



(d5f3 2 



p d5(3 2 
r dp 



2m~f 6/3 2 



27c 



dSfo _ p d5(3 2 
dr r dp 



r(l- p 2 ) 
mno6/3 2 



r(l - p 2 ) 



d 3 a 



(6.5) 



It is straightforward to perform these integrals numer- 
ically using the superfluid r-mode eigenfunctions de- 
scribed in Sec. and so determine the energy dissipa- 
tion rate caused by mutual friction. It is convenient to 
express this rate in terms of a mutual-friction damping 
time, t mf : 



1 1 fd£ 

Tmf 2£ \ dl 



1 

1~MF 



n 



yhrGpa 



(6.6) 



The characteristic mutual-friction damping time f MF , 
also defined in Eq. ( |6.6[ ) , is independent of angular veloc- 
ity and temperature (to lowest or der) . The f2 5 scaling of 
1/t mf follows directly from Eqs. (6.2) and flS.Ep. The Q 7 
scaling of d£ j dt in Eq. (|6.5[) follows in turn from Eq. (6.3) 



and the fact that 6w a scales as f2 3 to lowest order. We 
present in Figs. 6 and f^the numerically determined char- 
acteristic damping times t mf for the m — 2 superfluid 
r-modes. Fig. 6 illustrates how sensitively f MF depends 
on the entrainment parameter e. We see that f MF has a 
typical value of about 10 4 s, but that it becomes much 
smaller for a few narrow ranges of e. These spikes in 
the f MF (e) curves are caused by the "resonance" phe- 
nomenon that we discuss in Sec. |v|. The mutual friction 
damping time f MF becomes very small when the eigen- 
function 6(3 2 becomes large. This occurs at the points 
where the operator D has a vanishing eigenvalue. We see 
from Fig. || that the locations of these vanishing eigen- 
values coincide exactly with the location of the spikes 
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p s /10 14 g cm 3 

FIG. 7. Characteristic damping time ï M f due to superfluid 
mutual friction. This curve shows the dependence of f MF as 
a function of the superfluid transition density p a for e = 0.04. 



in the p s = 2.8 x 10 14 g/cm 3 curve in Fig. 6. We find 
that near these spikes the curve f MF (e) oc (e— e c ) 2 . This 
quadratic dependence is exactly what is expected given 
that Sp2 oc 'i/'min/Amin near these spikes. Fig. 6 also illus- 
trates that the exact location of these spikes depends on 
the value of the superfluid transition density p s . Fig. |^ 
illustrates that f MF changes smoothly as p s is varied: the 
locations of these spikes move smoothly to larger values 
of e as p s increases. Fig. 6 also illustrates that the prob- 
ability of having a small f MF does not depend strongly 
on p s . For example, we find that about 1% of the val- 
ues of e in the acceptable range, 0.02 < e < 0.06, have 
t M f(c) < 5 s, and this percentage is relatively insensitive 
to p s . 

The other important internal fluid dissipation mecha- 
nism in superfluid neutron stars is expected to be regu- 
lar shear viscosity. In the superfluid core this viscosity 
is due to clcctron-clectron scattering, while in the sur- 
rounding ordinary-fluid envelope the Standard neutron- 
neutron scattering dominates. The rate at which energy 
is dissipated by shear viscosity is given by: 



d£ 
~dï 



2 V Sa* ab Sa ab d^ 



2i 1 n a 5v b 5a* ab d 2 x. 



(6.7) 



The volume integral in Eq. (|6.7| ) is to be evaluated within 
the superfluid core, and within the ordinary-fluid enve- 
lope, but not over the boundary surface between the two. 
The surface integrals are to be evaluated over the interior 
(superfluid side) and exterior (ordinary-fluid side) of the 
superfluid boundary respectively. The tensor 5a a b that 
appears in these integrals is the shear of the electron ve- 
locity of Eq. (2.4) within the superfluid core, and the 
ordinary-fluid velocity öv a in the envelope. In the small 
angular velocity expansion the electron velocity öv® is 
the same as Sv a , to lowest order. Thus to lowest order 



the energy dissipation rate due to shear viscosity is given 
by 



d£ 

—— I = — ar 
dt I to 



,tt(to + l) 3 



(to- l)(2m+ l)!i?j^ 2 x 



Ra / \ 2m 



(2TO+1) J 7] I dr 



R, 



2m 



-(r l .-ri )R a \^j j, (6.8) 

where rj s and r\ are the limits of the viscosity taken from 
the superfluid and the ordinary-fluid side of the boundary 
respectively. 

In the superfluid core of the neutron star, r < R s , the 
appropriate viscosity to use in Eq. ( 3.5 ) is due to electron- 
electron scattering. This electron-electron scattering vis- 
cosity is given approximately by the expression J37| 



t? = 6.0x 10 6 (^)' 



(6.9) 



where all quantities are given in cgs units. Similarly, in 
the ordinary-fluid envelope, neutron-neutron scattering 
viscosity dominates at the densities where most of the 
dissipation occurs. This neutron-neutron scatte ring vis- 
cosity is given approximately by the expression [|37f 



7] = 347- 



3 9/4 
y2 ' 



(6.10) 



Using these expressions for the viscosity it is straight- 
forward to evaluate the energy dissipation rate for shear 
viscosity using Eq. (3.8). As in the case of mutual fric- 



tion, it is convenient to express the result as a viscous 
time-scale: 



2S\~dt 



1 /10 9 K 



T 



(6.11) 



27 1 n a Sv b e 5a* ab d 2 x P 



The chara cteristic viscous time-scale f v , also defined in 
Eq. ( 6.11 ), is independent of the angular velocity and 
the temperature of the neutron star. We find that f v 
1.01 x 10 8 s for our superfluid neutron-star model with 
= 2.8 x 10 14 g/cm 3 . This value is somewhat shorter 



than that obtained for hot neutron stars, 2.52 x 10 s, 
by Lindblom, Owen, and Morsink fï], and for cold neu- 
tron stars, 2.25 x 10 8 s, by Andersson, Kokkotas, and 
Schutz §. 

Gravitational radiation is the final form of dissipation 
expected to have a significant influence on the r-modcs 
of superfluid neutron stars. Since the superfluid r-modcs 
are identical to their ordinary-fluid counterparts to lowest 
order in the angular velocity, the gravitational radiation 
coupling is the same to this order. The characteristic 
gravitational radiation time-scale f GR for the to = 2 r- 
mode, 



1 



1 

28 



d£ 
17 



ï 



n 



(6.12) 
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therefore has the same value, f GR = 3.26 s, as in the 
ordinary-fluid case || . 

The effects of mutual friction, shear viscosity, and grav- 
itational radiation act together simultaneously to influ- 
ence the evolution of the superfluid r-modes. Their com- 
bined effects on the evolution of the energy of the mode 
are conveniently described by the overall dissipative time- 
scale, t: 



1 




d£ 
~dï 



(6.13) 



Since the energy £ is a real functional, the quantity 1/r 
is in fact just the imaginary part of the frequency of the 
mode. The sign of r therefore determines whether the 
mode is stable. If r > then dissipation decreases the 
energy of the mode and it is stable, while if r < then 
dissipation causes the energy (and hence the mode itself) 
to increase exponentially. We can explore the conditions 
under which mutual friction and viscosity are effective in 
suppressing the gravitational radiation driven instability 
by giving the explicit temperature and angul ar-vc locity 
dependence of r using Eqs. (3.6), (6.11), and (6.12): 



1 



1 



r(Ü,T) 



1~GR \ V TïGpo 

1 f Q 



10 9 K 
T 



(6.14) 



Gravitational radiation tends to drive the r-modes un- 
stable while mutual friction and shear viscosity tend to 
stabilize these modes. From Eq. (6.14) we see that r > 0, 
and mutual friction will completely suppress the gravita- 
tional radiation driven instability whenever 



< ï GR . 



(6.15) 



Since neutron stars have angular velocities that are lim- 
ited by f2 < ^^/irGpo, we see that mutual friction will 
suppress the gravitational radiation instability for all an- 
gular velocities whenever f MF < 1.5r GH ss 4.89 s. From 
Fig. 6 we see that this may occur, but only if the en- 
trainment parameter e of neutron-star matter is limited 
to a very narrow range. Only about 1% of the values of e 
in the expected range, 0.02 < e < 0.06, have sufficiently 
short f MF to supp ress the instability completely. 

From Eq. (|6.15| ) we see that mutual friction will always 
suppress the gravitational radiation driven instability in 
the r-modes of neutron stars with sufficiently small an- 
gular velocities. However, this suppression is not physi- 
cally relevant if f MF i s too large. If the angular velocities 
needed in Eq. ( 6.15 ) are sufficiently small, then shear 



viscosity will also play a significant role in suppressing 
the r-mode instability in these stars. In particular shear 
viscosity alone will suppress the r-mode instability if 
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FIG. 8. Critical angular velocities for superfluid neu- 
tron stars having a range of characteristic mutual friction 
time-scales t M f- 
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(6.16) 



Mutual friction then is capable of further suppressing the 
r-mode instability in more rapidly rotating neutron stars 
only if 



Tgi- 



1/6 / rp x 1/3 

\Wk) 



57Js (ïcfe) 



1/3 



(6.17) 



We see that mutual friction will be primarily respon- 
sible for the suppression of the r-mode instability in 
some (sufficiently warm) superfluid neutron stars only 
if t mf < 57.7 s. Such small values of f MF occur in 
only about 3% of the values of e in the expected range: 
0.02 < e < 0.06. 

It is also informative to examine the critical angular 
velocities, fi c , that mark the dividing line between stable 
and unstable neutron stars: t(SI c , T) — 0. Stars rotating 
more rapidly than Q c are unstable, while thosc rotating 
more slowly are stable. Fig. || illustrates the temperature 
dependence of Q c for a range of possible mutual friction 
time-scales. From Fig. || we see that shear viscosity com- 
pletely suppresses the gravitational radiation instability 
in all neutron stars cooler than about 10 6 K. We also 
see that the mutual friction time-scale must be shorter 
than 10 or 15 s in order for mutual friction to play a 
significant role in suppressing the gravitational radiation 
instability in neutron stars with temperatures that are 
typical of low mass x-ray binaries (about 10 8 K). Since 
only about 2% of the expected range of e have time-scales 
this short, it appears unlikely that mutual friction is act- 
ing to suppress the gravitational radiation instability in 
these stars. 
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APPENDIX: PULSATION EQUATIONS IN 
SPHERICAL COORDINATES 



In order to solve the superfluid pulsation equations nu- 
merically, it is necessary to express them in some partic- 
ular coordinate representation. We find it useful to work 
in spherical coordinate: r, p = cosö, and cp. The oper - 
ators D, E, F, D, E and F that appear in Eqs. ( 4.1S| ) 
and ( 1.19] ) have the following representations in spherical 
coordinates, 
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2p döU 2 
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r dr 
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We note that the functions p 22 and h 22 that appear 
on the right sidcs of Eqs. ( [A3] ) and (A6) are parts 
of the second-order expansions of these quantities. In 
general the second-order density function has the form: 
p 2 {r,p) = p 20 (r) + p 22 {r)P 2 {p), where p 2 o{r) and p 22 (r) 
are functions of r and P 2 (p) is the l = 2 Legendre poly- 
nomial. A similar expression is satisfied by h 2 . These 
functions may be determined using the techniques de- 
scribed in Lindblom, Mendell, and Owen p3|. 
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